*************************************************
* Purpose -- Purpose -- TWFE analysis using 
* CBP data (as provided in JNR's replication
* files)
************************************************
clear all
* RUN 00_path_master.do FIRST TO GET FILEPATHS


************************************************
* Loading different samples
************************************************
* 1. CZ-state
frame create twfe1
frame change twfe1

use "$clean/cbp/cz_state_sample.dta", clear

* creating frames for sub-periods

frame copy twfe1 twfe1_00

frame twfe1_00: keep if inrange(year,2000,2016)


* 2. MSCZ

frame create mscz
frame change mscz

use "$clean/cbp/mscz_sample.dta", clear

* creating frames for sub-periods
frame copy mscz mscz_00

frame mscz_00: keep if inrange(year,2000,2016)


* 3. BCP

frame create bcp
frame change bcp

use "$clean/cbp/bcp_sample.dta", clear

* creating frames for sub-periods
frame copy bcp bcp_00

frame bcp_00: keep if inrange(year,2000,2016)


************************************************
* Creating different frames for distributed lag regressions
************************************************
frame copy mscz mscz_dl
frame change mscz_dl

* getting extra data so that no observations are dropped due to missing leads/lags

forval i = 2017/2020 {
	insobs 1
	replace year = `i' if _n == _N
}

forval i = 1985/1989 {
	insobs 1
	replace year = `i' if _n == _N
}

fillin pair_id_czstate year
drop if mi(pair_id_czstate)
bysort pair_id_czstate (statefips): replace statefips = statefips[1] if mi(statefips)
bysort pair_id_czstate (pair_id_num): replace pair_id_num = pair_id_num[1] if mi(pair_id_num)
bysort pair_id_czstate (czone): replace czone = czone[1] if mi(czone)

* Merge with minimum wage data from VZ
merge m:1 statefips year using "$raw/other/mw_state_annual_vz", keep(3) keepusing(max_mw) nogen
replace mw = max_mw if inrange(year,2017,2020) | inrange(year,1985,1989)
replace lmw = log(mw)
drop max_mw
sort pair_id_czstate year

frame copy bcp bcp_dl
frame change bcp_dl

forval i = 2017/2020 {
	insobs 1
	replace year = `i' if _n == _N
}

forval i = 1985/1989 {
	insobs 1
	replace year = `i' if _n == _N
}

fillin pair_id_county year
drop if mi(pair_id_county)
bysort pair_id_county (statefips): replace statefips = statefips[1] if mi(statefips)
bysort pair_id_county (county): replace county = county[1] if mi(county)
bysort pair_id_county (pair_id_num): replace pair_id_num = pair_id_num[1] if mi(pair_id_num)

* Merge with minimum wage data from VZ
* For TWFE with CBP, we use JNR's data for their sample period (to replicate their results), and only supplement with Vaghul and 
* Zipperer data downloaded by us for out of sample years
merge m:1 statefips year using "$raw/other/mw_state_annual_vz", keep(3) keepusing(max_mw) nogen
replace mw = max_mw if inrange(year,2017,2020) | inrange(year,1985,1989)
replace lmw = log(mw)
drop max_mw
sort pair_id_county year

************************************************
* TWFE regressions
************************************************
foreach frame in  twfe1_00 twfe1 {
	frame `frame' {
		eststo `frame'_earn_w: reghdfe learn75 lmw learn_other lpop [aw=pop], a(czstate year) cluster(statefip)
		eststo `frame'_emp_w: reghdfe lemp75 lmw lemp_other lpop [aw=pop], a(czstate year) cluster(statefip)
		eststo `frame'_owe_w: ivreghdfe lemp75 (learn75 = lmw) learn_other lemp_other lpop [aw=pop], a(czstate year) cluster(statefip)
		estadd local fe " "
		estadd local czp_fe " "
		estadd local bcp_fe " "
		estadd local gap " "
		estadd local samples " "
		estadd local cz_state_sample "Y"
		estadd local mscz_sample " "
		estadd local bcp_sample " "

		eststo `frame'_earn_uw: reghdfe learn75 lmw learn_other lpop, a(czstate year) cluster(statefip)
		eststo `frame'_emp_uw: reghdfe lemp75 lmw lemp_other lpop, a(czstate year) cluster(statefip)
		eststo `frame'_owe_uw: ivreghdfe lemp75 (learn75 = lmw) learn_other lemp_other lpop, a(czstate year) cluster(statefip)
	}
}

foreach frame in  bcp_00 bcp {
	frame `frame' {	
		eststo `frame'_earn_w: reghdfe learn75 lmw learn_other lpop [aw=pop], a(pair_id_county pair_id_num##year) cluster(statefip)
		eststo `frame'_emp_w: reghdfe lemp75 lmw lemp_other lpop [aw=pop], a(pair_id_county pair_id_num##year) cluster(statefip)
		eststo `frame'_owe_w: ivreghdfe lemp75 (learn75 = lmw) learn_other lemp_other lpop [aw=pop], a(pair_id_county pair_id_num##year) cluster(statefip bordersegment)
		estadd local fe " "
		estadd local czp_fe " "
		estadd local bcp_fe "Y"
		estadd local gap " "
		estadd local samples " "
		estadd local cz_state_sample " "
		estadd local mscz_sample " "
		estadd local bcp_sample "Y"	
		
		eststo `frame'_earn_uw: reghdfe learn75 lmw learn_other lpop, a(pair_id_county pair_id_num##year) cluster(statefip)
		eststo `frame'_emp_uw: reghdfe lemp75 lmw lemp_other lpop, a(pair_id_county pair_id_num##year) cluster(statefip)
		eststo `frame'_owe_uw: ivreghdfe lemp75 (learn75 = lmw) learn_other lemp_other lpop, a(pair_id_county pair_id_num##year) cluster(statefip bordersegment)
	}
}

foreach frame in  mscz_00 mscz {
	frame `frame' {
		eststo `frame'_earn_w: reghdfe learn75 lmw learn_other lpop [aw=pop], a(pair_id_czstate pair_id_num##year) cluster(statefip)
		eststo `frame'_emp_w: reghdfe lemp75 lmw lemp_other lpop [aw=pop], a(pair_id_czstate pair_id_num##year) cluster(statefip)
		eststo `frame'_owe_w: ivreghdfe lemp75 (learn75 = lmw) learn_other lemp_other lpop [aw=pop], a(pair_id_czstate pair_id_num##year) cluster(statefip)
		estadd local fe " "
		estadd local czp_fe "Y"
		estadd local bcp_fe " "
		estadd local gap " "
		estadd local samples " "
		estadd local cz_state_sample " "
		estadd local mscz_sample "Y"
		estadd local bcp_sample " "	
		
		eststo `frame'_earn_uw: reghdfe learn75 lmw learn_other lpop, a(pair_id_czstate pair_id_num##year) cluster(statefip)
		eststo `frame'_emp_uw: reghdfe lemp75 lmw lemp_other lpop, a(pair_id_czstate pair_id_num##year) cluster(statefip)
		eststo `frame'_owe_uw: ivreghdfe lemp75 (learn75 = lmw) learn_other lemp_other lpop, a(pair_id_czstate pair_id_num##year) cluster(statefip)
	}
}

* Output twfe table
esttab twfe1_earn_uw mscz_earn_uw  bcp_earn_uw twfe1_00_earn_uw mscz_00_earn_uw  bcp_00_earn_uw ///
using "$tables/Table 1.tex", replace booktabs b(3) se(3) nomtitles keep(lmw) ///
varlabels(lmw "Log wages") nonote noobs star(* 0.10 ** 0.05 *** 0.01) width(\hsize) ///
mgroups("1990-2016" "2000-2016", pattern(1 0 0 1 0 0) span ///
prefix(\multicolumn{@span}{c}{) suffix(}) erepeat(\cmidrule(lr){@span})) ///
title(MW effects on log average earnings and log employment in restaurants: TWFE estimates for different samples and specifications using CBP data\label{tab:twfecbp}) ///
postfoot( ) ///
posthead(\midrule \\ ///
\multicolumn{7}{l}{\textbf{A. Unweighted}} \\ & & & \\)

esttab twfe1_emp_uw mscz_emp_uw  bcp_emp_uw twfe1_00_emp_uw mscz_00_emp_uw  bcp_00_emp_uw  ///
using "$tables/Table 1.tex", append booktabs b(3) se(3) nomtitles keep(lmw) ///
varlabels(lmw "Log employment") nonote noobs nonum star(* 0.10 ** 0.05 *** 0.01) ///
width(\hsize) prehead( ) posthead( ) postfoot( )

esttab twfe1_owe_uw mscz_owe_uw  bcp_owe_uw twfe1_00_owe_uw mscz_00_owe_uw bcp_00_owe_uw using ///
"$tables/Table 1.tex", append booktabs b(3) se(3) nomtitles keep(learn75)  ///
varlabels(learn75 "OWE") nonote noobs nonum star(* 0.10 ** 0.05 *** 0.01) width(\hsize) ///
prehead( ) posthead( ) postfoot( ) ///

esttab twfe1_earn_w mscz_earn_w bcp_earn_w twfe1_00_earn_w mscz_00_earn_w bcp_00_earn_w  using ///
"$tables/Table 1.tex", append booktabs b(3) se(3) nomtitles keep(lmw) ///
varlabels(lmw "Log wages") nonote noobs nonum star(* 0.10 ** 0.05 *** 0.01) width(\hsize) ///
prehead(\multicolumn{6}{l}{\textbf{B. Weighted by population}} ///
\\ & & & \\) postfoot( ) posthead( )

esttab twfe1_emp_w mscz_emp_w bcp_emp_w twfe1_00_emp_w mscz_00_emp_w bcp_00_emp_w  using ///
"$tables/Table 1.tex", append booktabs b(3) se(3) nomtitles keep(lmw) ///
varlabels(lmw "Log employment") nonote noobs nonum star(* 0.10 ** 0.05 *** 0.01) ///
width(\hsize) prehead( ) posthead( ) postfoot( )

esttab twfe1_owe_w mscz_owe_w bcp_owe_w twfe1_00_owe_w mscz_00_owe_w bcp_00_owe_w  using ///
"$tables/Table 1.tex", append booktabs b(3) se(3) nomtitles keep(learn75) ///
varlabels(learn75 "OWE") nonote noobs nonum star(* 0.10 ** 0.05 *** 0.01) width(\hsize) ///
stats(N gap fe czp_fe bcp_fe gap samples cz_state_sample mscz_sample bcp_sample, ///
labels("\textbf{N}" " " "\textbf{Fixed effects}" "MSCZ pair-year" "BC pair-year" ///
 " " "\textbf{Sample}" "All CZ-state" "MSCZ pairs" "Border county pairs") ///
fmt(%9.0fc)) /// 
prehead( ) posthead( ) ///
postfoot( ///
\hline\hline ///
\end{tabular*} ///
{\centering ///
\caption*{\begin{footnotesize} ///
\textit{Notes:} Estimated using specification described in equation \ref{eq:TWFE_logMW} and County Business Patterns data. Columns (1)-(3) use the years 1990-2016, and (4)-(6) use 2000-2016. Columns (1) and (4) use the sample of all commuting-zones-by-state, and include CZ-state and year fixed effects. Columns (2) and (5) use multi-state commuting-zone-by-state pairs sample and include fixed effects for CZ-state (repeated if a CZ-state is part of multiple pairs) and pair-year fixed effects. Columns (3) and (6) similarly use the sample of border county pairs and include fixed effects for border counties (each time a county is in a pair) and pair-year fixed effects. Outcomes include log weekly earnings, log employment, and own wage elasticity (OWE). OWE is the IV estimate from regressing log employment on log earnings instrumented by log minimum wage. $\dagger$ indicates OWE estimates with a weak first-stage. Estimates in Panel A are unweighted, while estimates in Panel B are weighted using the working-age population. All specifications include log of earnings/log of employment outside the restaurant sector, and working-age population as control variables. Standard errors, reported in parentheses, are clustered at the state level. Stars indicate statistical significance as follows: \\ ///
\sym{*} \(p<0.10\), \sym{**} \(p<0.05\), \sym{***} \(p<0.01\) ///
\end{footnotesize}}} ///
\end{table} ///
) ///
substitute("\_" "_")


***************************************************
* TWFE DL (recreate JNR table 7)
***************************************************
eststo clear
local f = 1
foreach frame in bcp_dl mscz_dl {
	if `f' == 1 {
		local a = "pair_id_county"
	}
	
	if `f' == 2 {
		local a = "pair_id_czstate"
	}
	
	frame change `frame'
	
	gen lmw_f_3_1 = F3.lmw - F.lmw
	gen lmw_f_1_0 = F.lmw - lmw
	gen lmw_f_4_2 = F4.lmw - F2.lmw
	gen lmw_f_2_0 = F2.lmw - lmw
	
	reghdfe lemp75 lmw_f_3_1 lmw_f_1_0 lmw lemp_other lpop if year <= 2013, ///
	a(`a' pair_id_num##year) cluster(statefips)
	gen flag1 = e(sample)
	
	eststo cbp_`frame'_3_1: nlcom (t_3: _b[lmw_f_3_1]) (t_1: _b[lmw_f_1_0]) ///
	(t0: _b[lmw]) (trend: _b[lmw_f_1_0] - _b[lmw_f_3_1]), post
	
	reghdfe lemp75 lmw_f_4_2 lmw_f_2_0 lmw lemp_other lpop, ///
	a(`a' pair_id_num##year) cluster(statefips)
	gen flag2 = e(sample)
	
	eststo cbp_`frame'_4_2: nlcom (t_4: _b[lmw_f_4_2]) (t_2: _b[lmw_f_2_0]) ///
	(t0: _b[lmw]) (trend: _b[lmw_f_2_0] - _b[lmw_f_4_2]), post
	
	local ++f
}

esttab cbp_mscz_dl_3_1 cbp_bcp_dl_3_1 using "$tables/Table A1.tex", replace ///
booktabs b(3) se(3) mtitles("MSCZ pairs" "BC pairs") varlabels(t_3 "t: -3" ///
t_1 "t: -1" t0 "t: 0" trend "Trend (t[-1] - t[-3])") star(* 0.1 ** 0.05 *** 0.01) width(\hsize) ///
title(Modification of JNR table 7 \label{tab:jnr_tab_7}) noobs nonotes ///
postfoot( ) ///
posthead(\midrule \\ ///
\multicolumn{3}{l}{\textbf{A. Trend between t: -3 and t: -1}} \\ & & & \\) ///
substitute("\_" "_")

esttab cbp_mscz_dl_4_2 cbp_bcp_dl_4_2 using "$tables/Table A1.tex", append ///
booktabs b(3) se(3) varlabels(t_4 "t: -4" t_2 "t: -2" t0 "t: 0" trend "Trend (t[-2] - t[-4])") ///
star(* 0.1 ** 0.05 *** 0.01) width(\hsize) noobs nonum nomtitles ///
prehead(\multicolumn{3}{l}{\textbf{B. Trend between t: -4 and t: -2}} \\ & & & \\) ///
posthead( ) ///
postfoot( ///
\hline\hline ///
\end{tabular*} ///
{\centering ///
\caption*{\begin{footnotesize} ///
\textit{Notes}: Panel A estimated using the specification in Equation 7 in JNR. Panel B is a modification of the equation with the third and first lag replaced with the fourth and second, respectively. Both columns use data from County Business Patterns. JNR do not augment their data with out-of-sample minimum wage observations, so some years get dropped because of missing leads of log minimum wage. Thus, their results contain observations up till 2013. To ensure replicability of Panel A, we also restrict to observations until 2013. For Column B, as we have an additional lead, we augment with minimum wage data from 2017, so that we can still use data up to 2013. The distributed lag figures in Figure \ref{fig:dl_uw}, on the other hand, use our preferred method of using all available data until 2016. In this table, Column 1 uses multi-state commuting-zone-by-state pairs sample and column 2 uses the sample of border county pairs. Standard errors, reported in parentheses, are clustered at the state level. Stars indicate statistical significance as follows: sym{*} \(p<0.10\), \sym{**} \(p<0.05\), \sym{***} \(p<0.01\) ///
\end{footnotesize}}} ///
\end{table} ///
) ///
substitute("\_" "_")



***************************************************
* TWFE DL figure
***************************************************

frame change mscz_dl

reghdfe lemp75 L(1/5).lmw lmw F(1/3).lmw lemp_other lpop, a(pair_id_czstate pair_id_num##year) cluster(statefips)

nlcom (t_4: -(_b[F1.lmw] + _b[F2.lmw] + _b[F3.lmw])) ///
(t_3: -(_b[F1.lmw] + _b[F2.lmw])) (t_2: -_b[F1.lmw]) (t_1: 0) ///
(t0: _b[lmw]) (t1: _b[lmw] + _b[L1.lmw]) ///
(t2: _b[lmw] + _b[L1.lmw] + _b[L2.lmw]) ///
(t3: _b[lmw] + _b[L1.lmw] + _b[L2.lmw] + _b[L3.lmw]) ///
(t4: _b[lmw] + _b[L1.lmw] + _b[L2.lmw] + _b[L3.lmw] + _b[L4.lmw]) ///
(t5: _b[lmw] + _b[L1.lmw] + _b[L2.lmw] + _b[L3.lmw] + _b[L4.lmw] + _b[L5.lmw])

matrix step1 = r(table)[1,1..10]
matrix step2 = r(table)[5..6,1..10]
matrix lemp = (step1 \ step2)
matrix lemp = lemp'

frame create mscz_dl_lemp
frame change mscz_dl_lemp
svmat lemp, names(col)
foreach v in ll ul {
	replace `v'=0 if mi(`v')
}
gen eventtime = _n - 5


frame mscz_dl_lemp: twoway (line b eventtime, lcolor(blue*0.5) lpattern(solid)) ///
(rcap ll ul eventtime, lcolor(blue*0.5) lpattern(dash)), ///
xtitle(Years around MW increase event) xlabel(-4 "{&le} -4" -3(1)4 5 "{&ge} 5") ///
ytitle("Cumulative response of log employment") legend(off) ///
yline(0) ylabel(-0.8(0.2)0.4) title(B. Sample - CBP) ///
saving("$figures/mscz_cbp_dl", replace)

***************************************************
* TWFE rolling regression figures
***************************************************

frame change mscz

local j = 2016 
	forvalues i = 1990/2005 {
	    
		reghdfe lemp75 lmw lemp_other lpop [aw=pop]  if inrange(year,`i',`j'), a(pair_id_czstate pair_id_num##year) cluster(statefip)
		matrix w_cbp_`i' = (r(table)[1,1], r(table)[5,1],r(table)[6,1])
		
		if `i' == 1990 {
			matrix w_estimates`j' = (w_cbp_`i')
		}
		
		else {
			matrix w_estimates`j'_`i' = (w_cbp_`i')
			matrix w_estimates`j' = (w_estimates`j' \ w_estimates`j'_`i')
		}
		
		reghdfe lemp75 lmw lemp_other lpop [aw=lpop]  if inrange(year,`i',`j'), a(pair_id_czstate pair_id_num##year) cluster(statefip)
		matrix w_lpop_cbp_`i' = (r(table)[1,1], r(table)[5,1],r(table)[6,1])
		
		if `i' == 1990 {
			matrix w_lpop_estimates`j' = (w_lpop_cbp_`i')
		}
		
		else {
			matrix w_lpop_estimates`j'_`i' = (w_lpop_cbp_`i')
			matrix w_lpop_estimates`j' = (w_lpop_estimates`j' \ w_lpop_estimates`j'_`i')
		}
		
		
		reghdfe lemp75 lmw lemp_other lpop if inrange(year,`i',`j'), a(pair_id_czstate pair_id_num##year) cluster(statefip)
		matrix uw_cbp_`i' = (r(table)[1,1], r(table)[5,1],r(table)[6,1])
		
		if `i' == 1990 {
			matrix uw_estimates`j' = ( uw_cbp_`i')
		}
		
		else {
			matrix uw_estimates`j'_`i' = (uw_cbp_`i')
			matrix uw_estimates`j' = (uw_estimates`j' \ uw_estimates`j'_`i')
		}
		
	}	

foreach mat in w_estimates2016 w_lpop_estimates2016 uw_estimates2016 {
	matrix colnames `mat' = b_cbp ll_cbp ul_cbp	
}

foreach i in w w_lpop uw {
	frame create estimates_`i'
	frame change estimates_`i'
	svmat `i'_estimates2016, names(col)
	gen year = _n + 1989
	save "$intermediate/rreg_cbp_estimates_`i'.dta", replace
}



